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Ab-Initio Theory of Superconductivity in a Magnetic Field II. : Numerical solution. 
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We numerically investigate the Spin Density Functional theory for superconductors (SpinSCDFT) 
and the approximated exchange-correlation functional, derived and presented in the preceding paper 
I. As a test system we employ a free electron gas featuring an exchange-splitting, a phononic pairing 
field and a Coulomb repulsion. SpinSCDFT results are compared with Sarma, the Bardeen Cooper 
and Schrieffer theory and with an Eliashberg type of approach. We find that the spectrum of the 
superconducting Kohn-Sham SpinSCDFT system is not in agreement with the true quasi particle 
structure. Therefore, starting from the Dyson equation, we derive a scheme that allows to compute 
the many body excitations of the superconductor and represents the extension to superconductivity 
of the Go Wo method in band structure theory. This superconducting Go Wo method vastly improves 
the predicted spectra. 


I. INTRODUCTION 

Interaction between the magnetic (M) and supercon¬ 
ducting (SC) order leads to complex and fascinating phe¬ 
nomena. Apart from the Meissner effect as the most ap¬ 
parent aspect of this interaction on macroscopic length 
scales, for singlet superconductors, the ferromagnetic 
parallel spin alignment competes with spin anti-parallel 
Cooper pair formation. While for triplet superconduc¬ 
tors such as UGe 2 a ferromagnetic (F) order is possible 
even in a bulk geometry 1 , F/SC interfaces or SC sur¬ 
faces in an external magnetic field allow to study the 
microscopic competition of a large spin-splitting also for 
singlet superconductors. This may lead to spatial in¬ 
homogeneities of the SC order parameter, such as the 
phase predicted by Fulde and Ferrell and Larkin and 
Ovchinnikov 2,3 . Furthermore, the spin valve behavior 
of complex F/SC structures 4 6 may provide opportuni¬ 
ties for novel devices making use of the unique electronic 
configuration that appears due to the vicinity of these 
two competing phases (see Ref. 7 for a review). 

These effects are addressed in the theoretical literature 
so far mostly within model or semi-empirical calculations 
due to the lack of a complete and efficient ab-initio the¬ 
ory. This leaves the prediction of essential material de¬ 
pendent properties as critical temperature and excita¬ 
tion gap in the presence of a magnetic field out of reach. 
The Spin Density Functional theory for superconductors 
(SpinSCDFT) approach presented by Ref. 8 (hereafter 
referred to as I) may fill this gap, as the theory has the 
computational convenience of a Kohn-Sham density func¬ 
tional framework and allows to calculation of material de¬ 
pendent SC parameters from the crystal structure. The 
SpinSCDFT is in principle exact, but relies on the ap¬ 
proximation of the exchange-correlation (xc) potential. 
A first approach to derive such an xc potential relies, in 
turn, on the Sham-Schluter equation 9 for a SC 10 and is 
presented in I. 

In this work, we present numerical results for Spin¬ 
SCDFT, aiming to achieve a deeper understanding of 
this theoretical framework and to characterize and vali¬ 
date the xc potential, as derived in I. In particular we 


will investigate the properties of the Kohn-Sham pair¬ 
ing function, that is a key object in SCDFT 11,12 . The 
test system we adopt for this analysis is a spin-splitted 
free electron gas with a phononic and Coulomb coupling. 
Details of the model will be presented in Sec. II. One 
advantage of this simplified model with a homogeneous 
exchange splitting is its similarity to the starting point 
of Ref. 13 and 14 for their discussion of the Eliashberg 
equations and BCS theory, respectively. As compared 
to Ref. 13, we use a different notation (compare I) and 
take a more general route which reduces to the earlier re¬ 
sults in the case that the magnetic field homogeneously 
splits the electronic states. We will compute the temper¬ 
ature vs exchange splitting diagram of the model using, 
apart from SpinSCDFT, the BCS theory and the Eliash¬ 
berg equations. Then, in Sec. V we will compare our 
SpinSCDFT results with the BCS approach (reviewed 
in Sec. Ill) and with the reference Eliashberg method 
(Sec. IV). 

The SpinSCDFT Kohn-Sham system proves to give 
qualitatively correct results for the J — T diagram. How¬ 
ever, we find in Sec. V that it does not show a physi¬ 
cal excitation spectrum. A similar problem is very well 
known in conventional DFT, and is usually called the 
band gap problem. Since the excitation gap is a very im¬ 
portant property of superconductors, it is important to 
devolve methods to compute it. Therefore, the last part 
of this work will be devoted to describe an extension of 
the GoWo method to our superconducting system and 
show that it entirely solves the problem, similar to its 
normal state counterpart 15 . 


II. A TEST SYSTEM 

The model system which we will use to investigate the 
SpinSCDFT formalism is based on a non interacting elec¬ 
tron gas under the influence of an homogeneous magnetic 
field B 0 . The energy of its electronic states relative 
to the Fermi energy Ef (k = k,n where here n is a band 
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Figure 1. (color online) Model c?F{uj) function used in this 
work (full red line) as compared to that of MgB 2 (green dotted 
line). 



Figure 2. (color online) Screened Coulomb interaction func¬ 
tion C stat (e,e / ) as given by the model expression in Eq. (3). 
The chosen model parameters are Ef = 4.78 Ha and k^ F — 
(0.005) 2 Ha. 


index and we use the notation —k = —fe, n), reads 

£ka = Ifc 2 - E f - sign(c7)/u B Bo • (1) 

The Fermi energy is defined by integrating the density of 
states (DOS) up to Ef to have N e electrons in the sys¬ 
tem. We set the density to N e /Q uc = 1 a^ 3 (ao is the 
Bohr radius and Q uc the unit cell volume) which leads 
to a relatively large Ef of 4.78 Ha ignoring the small 
imbalance in up and down spin occupations. We also 
define a center of energy between spin splitted states 
s(k) = \(ek^ + e -ki) = \k 2 — Ef and the splitting 
J(k) = \{ek\ — £-ki) = — EbB 0 . This will prove use¬ 
ful since, as seen in I, many SpinSCDFT entities depend 
on k only via these two parameters 6 and J. 

Superconductivity is induced in this model by an 
electron-phonon like attractive interaction, expressed by 
the Gaussian Eliashberg function 16 : 

o—, x . UJ 1 1 / uj — up \ 2 

c?F(u) = A—-—e “w ' . (2) 

2 ^w v 


This model depends on three parameters: A the electron- 
phonon coupling constant 16 ; ujq the center of mass of 
the phonon spectrum, and uo w the width of the opti¬ 
cal branch. In the calculations we fix these numbers to 
cjo = 2.2mHa ujw = 0.5mHa and A = 0.7 which lead 
to coupling properties that are loosely similar to those 
of MgB 2 17 . The resulting spectrum is plotted in Fig. 1 
and compared with a recomputed one of MgB 2 . In Spin¬ 
SCDFT one can consider a general Coulomb coupling 
on the same footing as the phonon interaction 11 . Here, 
we use a simple Thomas-Fermi based model that was 
used before in SCDFT 11,18 . In this model the screened 
Coulomb matrix elements between a state of energy 6 
and one of energy s' is given by 


-Ve • • Ef) 

£+£ / +2-E 1 / + 2\/(£+.E/)(£ / +-5'/) + 
£+£'+2 Ef - 2 y/{e+Ef){e'+Ef) + \k% v ) 


The screening parameter is chosen to be k^ F = 
(0.005) 2 Ha. With this parameter, the shape of the model 
C' stat (e,£ / ) is shown in Fig. 2 . All properties of the test 
system depend on the Bloch vector k and the band index 
n only via the single particle energy £&. For brevity we 
use the notation e = (e, J) f de = f ds f d J. Further 
let 6(e — e') = 5 [s — s')b{J — J'), then we may cast a Bril- 
louin zone integral into the isotropic formulation with the 
double DOS 


eM = £5 (E -£tty±t)S(j- £ -ii^i±). (4) 

k 

This quantities describes the number of states on equal 
center of energy 5 and splitting J surfaces. In our model 
p(e) the external field is homogeneous. This means the 
number of states on equal splitting surfaces has a delta 
distribution character that peaks at J 0 = — n F B 0 . In the 
remainder of the paper, J 0 replaces the J integrals almost 
everywhere so we simplify the notation using J 0 —>> J. 


III. THE BCS THEORY WITH AN EXCHANGE 
SPLITTING 


The J—T diagram of a BCS model with a homogeneous 
exchange splitting parameter J has been presented by 
Ref. 14. This approach, that we are going to review here, 
can only be used to obtain qualitative results. Still, it 
will be an important guideline in understanding the more 
involved Eliashberg and SpinSCDFT results of the next 
sections. In a BCS model 19 one replaces the interactions 
among single electrons with an effective one, keeping only 
the matrix elements that couple the states &, j" and — fc, 
The effective interaction is approximated with “a box” 
centered at the Fermi level (from — i? d to f2 d which is 
of the order of the Debye phonon frequency to mimic 
phononic type of pairing and with height —V). This 
leads to a fixed point equation for the mean field gap 
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A 14 

- fp(J + + /A 2 )). (5) 

p(0) is the DOS at the Fermi level and J is the splitting 
energy between up and down states. Apart from the 
solutions A of Eq. (5) there is also the trivial solution 
A = 0. We solve Eq. ((5)) numerically as a function of 
T and J 20 . The solutions A(T, J ) are presented in Fig. 3 
a). There, we normalize A to Z\ 0 , the solution for T —>> 0 
and J = 0. Similarly, we normalize the J to A 0 and T to 
T c0 , the critical temperature for J = 0. In this way we 
remove the explicit dependence on the parameters p(0)V 
and !? d . 

When one attempts to linearize Eq. (5), a peculiar be¬ 
havior is found in that the T C (J) curve bends inwards 14 . 
We solve the linearized Eq. (5) and show the resulting 
T C (J) as a green line in Fig. 3 a). As pointed out by 
the Refs. 14 and 23, unlike the original BCS model at 
J = 0, this equation leads to a J — T diagram in which 
the SC transition can be discontinuous in A, i.e. of first 
order. Below the temperature T/T c0 « 0.6 at point A 
and the dashed line in Fig. 3 b) no small A solution to 
the non-linear equation is can be found and the initial 
assumption of the linearization that an arbitrary small 
solutions exists is not valid. 

While we can find a non-vanishing solution A it may 
not correspond to the stable thermodynamic phase. In 
Fig. 3 b) we remove the non-vanishing solutions A , 
if the free energy favors the magnetic state. The re¬ 
sulting T — J diagram shows that for J larger than 
to the Chandrasekhar-Clogston limit at T = 0 21,22 of 
J c = A 0 / a/ 2 no SC solution is stable. 

Another interesting approach to describe SC in the 
presence of a magnetic field is presented by Powell et 
al. 23 who use a Hubbard model in connection with a ho¬ 
mogeneous exchange splitting. They treat the pairing 
part of the interactions among electrons in the system in 
the Hartree-Fock approximation, similar to BCS as de¬ 
scribed above and consequently arrive at a similar gap 
equation as compared to Eq. ((5)). The matrix elements 
of the KS system of SpinSCDFT within the spin decou¬ 
pling approximation will turn out to have a similar ana¬ 
lytic structure. 

Also, Ref. 23 discusses why the transition is of first or¬ 
der. They observe that for J < A and T = 0 the gap 
equation (5), and consequently A, is independent on J. 
Thus A(J) = A 0 for J < A 0 . At e < VJ 2 — A 2 on 
the other hand the Fermi functions at T = 0 in Eq. (5) 
are equal and thus cancel. Also for this type of solution 
A must be larger than J and for J > A only the triv¬ 
ial solution A = 0 can be found. At this point follows 
that A(T = OK, J) = Z\o0(Ao — J) and the transition is 
discontinuous T = 0. 

The above analysis will be crucial later, in Sec. IV and 
V, to guide the discussion of the more sophisticate ap¬ 
proaches, that feature a qualitatively similar behavior. 


In the next section we will discuss results of the Eliash- 
berg method (as derived in I, Sec. IV) when applied to 
our test system of Sec. II. 


IV. SOLUTIONS TO THE PHONON ONLY 
ELIASHBERG EQUATIONS 

We solve the Eliashberg Eqs. (1.130) to (1.132). The 
approximations used here, for the special case of homo¬ 
geneous exchange field, lead to equations similar to those 
derived by Vonsovsky et al. 13 . 

Similar to every equation that describes a sponta¬ 
neously broken symmetry, in addition to a possible fi¬ 
nite solution, the Eliashberg Eqs. (1.130) to (1.132) al¬ 
ways have the solution A^(J) = 0. Usually, this non-SC 
solution is not stable below T c in the sense that small 
symmetry breaking fields (that in the self consistent iter¬ 
ation scheme is equivalent to a small but non-zero start¬ 
ing guess) lead to the finite A^(J) solution via iteration 
of the Eliashberg Eqs. (1.130) to (1.132). Thus, we say 
that in this case the Z\®(J) = 0 solution has a zero basin 
of attraction; Only the starting value A ® init (J) = 0 leads 
to the final solution A^(J) = 0. Whenever J — 0, the 
Afi(J) = 0 solution has a zero basin of attraction below 
T 

C • 

From Eq. (1.113) we know that the complex A ® 
changes the poles of the Green function. We assume the 
term A£(cj n ) to be zero, for simplicity. Then from the 
analytic continuation to the real axis of Eq. (1.113), we 
see that the energy uj of such a pole satisfies the condi¬ 
tion L>j = sign (cr)Jk ± \J£k 2 + A(u;) 2 which is analogous 
to the usual Eliashberg equations (compare also Ref. 24) . 
At T = 0, the analytic continuation of A ® the real axis 
is purely real in the range of the Fermi energy and its 
value there defines the SC excitation gap 24 . Thus, the 
Matsubara component n = 0 of A® is related to the SC 
excitation gap of the quasi particle system. We choose 
this as a characteristic property that we investigate as a 
function of J and T. In the following we generate two 
J — T diagrams shown in Fig. 4. In a) we follow the SC 
solution, i.e. we take the converged A^(J) as input for 
the calculation at Z\^(J-f dJ), starting at J = 0 with dJ 
positive. This way we compute the diagram M from left to 
right" and test the stability of the A^(J) ^ 0 solution. 
In b), we take the converged Z\^(J) as input for the cal¬ 
culation at A^(J — dJ), starting at J = 0.5mHa. Thus, 
we generate the diagram "from right to left". Because for 
large J A^(J) is zero, we start from a small, symmetry 
breaking value at A%(J — dJ) instead of zero. This way 
we test the stability of the trivial A^(J) = 0 solution. 

Comparing a) and b) we see that the borders of sta¬ 
bility between the stability of A^(J) = 0 and A®(J) / 0 
do not agree. In fact, we find a region where both, the 
Al(J) = 0 and the A^(J) 0 solution have a finite 

basin of attraction; here the normal and the SC state 
are (meta) stable. The shape of the border of the region 
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Figure 3. (color online) BCS solutions for a spin splitted band structure 14 . In the panel a) we plot the solution A if we can find 
one, while in b) A is set to zero if the free energy favors the magnetic state. The green curve in a) shows the T C (J) behavior 
from the linearized equation which has a curious shape that bends inwards. Below the thin dashed line in b) at the label A at 
T/T c o ~ 0.6 no solution with small A exists and the transition is of first order. Label B at l/\/2 represents the Chandrasekhar- 
Clongston 21,22 limit. 


where J) = 0 is unstable resembles closely to the lin¬ 
ear BCS solution which we show in Fig. 4 a) as a green 
dashed line. 

We plot the Z\®(J) at T = 10K in Fig. 5 a) and 40K in 
Fig. 5 b) as a function of J on the vertical axis. The cor¬ 
responding equal temperature lines are blue in Fig. 4. We 
find that the shape is largely independent on the split¬ 
ting J and the temperature T except for a scale factor. 
Thus Z\^ =0 is sufficient to investigate the behavior of the 
theory. For low temperatures the down-scaling is much 
less pronounced and it is safe to say that the pairing is 
almost unaffected by the presence of a splitting up un¬ 
til the point where the SC phase is suppressed. For a 
high temperature, instead, the down-scaling is more pro¬ 
nounced and the transition becomes continuous above a 
certain temperature. 

V. RESULTS OF SPINSCDFT WITH THE 
GO-FUNCTIONAL 

In this Section we discuss the numerical solution of the 
SpinSCDFT gap equation (1.95) using the ire-potential 
derived in Sec. I.C.. We refer to this functional as the 
GO-functional. 

In Eq. (1.95) of I, we have derived the gap equation 
of SpinSCDFT using the GO-functional. This equations 
(1.95), in turn, is derived from the Sham-Schliiter equa¬ 
tion for a superconductor, written in I in the form 

Jdc , Sp[Al}(c,t , )Al(c') = 0 . (6) 

From the previous discussion in the Secs. Ill and IV, a 
continuous transition is to be expected for a small ex¬ 
change field intensity J as compared to the transition 
temperature. 

For the point of the continuous transition Eq. (6) can 
be linearized in A s s . Similar to I we use the notation with 


a breve to indicate linearized entities Sp = Sp [Z\g = 0]. 
Thus, in this case T C (J) can be computed from the con¬ 
dition that Sp = Sp + S™ + Sp has a singular eigenvalue 

detS/3 = 0 . (7) 

The corresponding shape of the solution Z\|/||Z\||| is the 
right eigenfunction to such a singular eigenvalue. 

Sp(t,t') is given in Eq. (1.80). To investigate the 
structure and properties of the SpinSCDFT xopotential 
is easier within the linearized form, since the matrix 
Sp(t, e') = *§§(e, e') + S^(e, e') + (e, e') is independent 

of the potential A s s . As discussed in detail in I, *§?(e, e') 
(5§(e, e')) corresponds to the Nambu (off) diagonal self¬ 
energy contribution. S™ is due to the v xc part of the 
Sham-Schliiter equation. In Sec. V A we present and dis¬ 
cuss the shape of the contributions *S^(e, e'), ^(e, e'), 

S^(e, e') and the T C (J) curve from the linearized xc- 
potential. 

Finally, the properties of the general non-linear gap 
equation, i.e. the J—T diagram of the solutions to Eq. (6) 
with and without the Coulomb repulsion will be pre¬ 
sented in Sec. VB. 


A. Linearized Sham-Schliiter Equation 

As discussed before, in the part of the J — T diagram 
for a relatively small applied field (i.e. low splitting J and 
high T) we expect a second order phase transition. This 
section deals with the corresponding continuous transi¬ 
tion. In Sec. VA1, will show the shape of 5^, and 

Sp . To determine the point of the transition according 
to Eq. (7), in Sec. VA2 we investigate the spectrum of 
Sp as a function of temperature and splitting and the 
corresponding solutions zA|/||zA|||. Then we will discuss 













5 




b) 


Figure 4. (color online) J — T diagram of the n = 0 component A E = o from the solution to the Eliashberg equations. We follow 
the SC solution A E (J) ^ 0 in a) or the non-SC solution A E (J) = 0 in b) and observe that we can find a region where both are 
(meta) stable. We show the full solution Z\ E (J) along the blue lines in Fig. 5. For comparison we show the linear BCS curve 
as a green dotted line in a). In b) we compare with Spin SCDFT results of Sec. V (green curve); the black curve is scaled on 
both axis by T c Elishberg /T c SpinSCDFT . 



Figure 5. (color online) Solutions to the Eliashberg equations 
Z\ E (J) for T = 40K a) and T = 10K b) along the blue lines 
of the left panel of Fig. 4., 


the shape of the T c ( J) curve in Sec. V A 3 from this linear 
approach. 


1. Temperature Dependence of Sp 


The three contributions to Sp(t, e') are (see Sec. I.C., 
Eqs. (1.91),(1.92),(1.93) and (1.94)) 

^(e,e') = (^(e) + ^(e))5(e-e') + 

+ ^ph/3( e > 0 + 0- (8) 

In this linear Sham-Schliiter form, S™(t) and Sp(t) mul¬ 
tiply Z\g(e) directly. They are shown for several T for 
J = O.OmHa and J = O.lmHa in Fig. 6 panel a) and b), 
respectively. Note the logarithmic center of energy scale 
e in all the plots in this section. The color scale (blue to 
red) indicates increasing temperatures. All terms have 
features only in the close vicinity to 5 = 0 and quickly 
decay to zero within a characteristic energy width of the 
phonon coupling. This energy scale is the analog of De¬ 
bey frequency u;o in Eq. (2). However the 5 dependence 


shown in Fig. 6 in the presence (panel b) and the ab¬ 
sence (panel a) of an exchange splitting is very different. 
In fact, in Fig. 6 a) where J = 0 both S™(t) and Sp (e) are 
positive and monotonously decreasing as a function of \s\. 
In presence of a J ^ 0 (Fig. 6 b), instead, they have the 
following complex temperature and energy dependence: 
For small T in the range \e\ < J, Sp (e) is negative and, in 

the limit T —>• 0, S™(t) and Sp (e) approach zero from op¬ 
posite sides. At \e\ ~ J both S^f(e) and Sp(t) vary very 
rapidly. This behavior is smoothed out with increasing T 
and at temperatures high enough with respect to J the 
non-splitted behavior is recovered. The temperature and 
J = OmHa and J = O.lmHa dependence of S^ h p(t, e 7 ) 

is shown in Fig. 6 c) and d), respectively. serves as 
a scale that other kernel contributions have to be com¬ 
pared with, so we choose a color scale that is relative to 
the maximum of S™, indicated on the right of every plot. 

For J = OmHa we note that the size of t f ) 

(Fig. 6 c) decays faster with temperature the one of the 
diagonal S™ and Sp (the the position of "white” in the 
color scale of Fig. 6 c) moves to the left with increasing 
temperatures). Furthermore, being both positive and di¬ 
agonal, Sp- and Sp have to be compared with the eigen¬ 
values of S^ h p. Sp- and Sp alone would result in a positive 
definite Sham-Schliiter matrix for J = OmHa (compare 
Fig. 6 a) at all temperatures so there is not non-trivial 
solution to Eq. (6). Thus, technically, the phase transi¬ 
tion from the SC to the non SC regime with the singular 
eigenvalue is induced by this relative reduction of S* h p 

as compared to S™ plus Sp. We will turn to a system¬ 
atic analysis of the eigenvalues of the linearized Sham- 
Schliiter matrix Sp in Sec. V A 2. 

The relative scale reduction is also found for the split- 
ted S^p. At e ~ 0, however, we stay much below the 
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(a) S^(t) (1.) and Sp (e) (r.) for J = OmHa 
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(c) *§p h/ 3 (e, e') for J = O.OmHa at (1. to r.) 

T = IK, 10K, 30K 






(b) S^(t) (1.) and Sp (e) (r.) for J — O.lmHa. 
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(d) Sp h p(z, e') for J — O.lmHa at (1. to r.) 
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(f) Scp(t, e') for J = 0.05mHa at (1. to r.) 

T = IK, 10K, 30K 


Figure 6. (color online) Contributions to the linearized Sham-Schliiter Eq. (6). In the top row we show the diagonal S™ and 
Sp that originate from the v xc and Nambu diagonal self-energy in the Sham-Schliiter equation, respectively. In the second 
(bottom) row, we show the contributions that originate from the Nambu off phonon (Coulomb) self-energy. The color scale 
of Sp(t,t') for negative values (decreasing: blue to white to green) is relative to max(S'^ 1 ) (white). Red to yellow to white 
indicates increasingly positive values. Note that Sp and S^ h p switch sign at e « 0 for J = O.lmHa at low T as compared to 
J — OmHa. 


scale of and exceed it only for higher temperatures. 
Moreover, the sign change of Sp is effectively reducing 
the diagonal repulsion. 

A purple line in Fig. 6 d) indicates the zero contour 
and shows that for very low T, S* h p(t, e') is positive for 
approximately the region where \e\ < J or \e'\ < J') and 
has a sharp negative spike at e = e' ^ J. Thus, as a 
curious fact, the phonon interaction is not "attractive" 
everywhere in this case. We show the shape of the static 
Coulomb part S*p in Fig. 6 row e) and f) for J = OmHa 
and J = 0.05mHa, respectively. Apart from the differ¬ 
ences in sign the overall behavior of the Coulomb term 
and phonon terms is roughly similar with significant de¬ 
viations in the fact that it does not change sign for a low 
temperature and exchange splitting, compare Fig. 6 d) 
with f). 

In summary we can say that we see relevant changes 
in the shape of the contributions to Sp for a finite ex¬ 
change splitting for the low temperature limits in the 
region |e:| < \J\ as compared to the spin degenerate case. 
At higher temperatures the splitting becomes less impor¬ 
tant. We point out that we know from the earlier discus¬ 
sion that this is the region, where we expect the lineariza¬ 


tion to be unjustified. From the form of the Bogo linbov 
eigenvalues E° = sign (cr)J + sign(<a )^£ 2 + |Z\| | 2 we ex¬ 
pect that, whenever A s s is larger than J, will see a behav¬ 
ior more similar to the case J — 0 . The reason is that, 
then, only the a = + branch has positive excitation en¬ 
ergies E+ > 0 , meaning that the ground state does not 
correspond to some of the excitations 7 ^ being occupied 
(see the discussion in I.III.A.2.dand by Ref. 14). 


2. Critical Temperatures and the Shape of Al 

Since we compute the critical temperature from Eq. 7, 
i.e. the occurrence of a singular eigenvalue of Sp, in this 
section, we will investigate the full spectrum as a func¬ 
tion of T and J. The KS potential A s s is proportional to 
the right eigenvector of Sp that is associated to a singu¬ 
lar eigenvalue. Thus, all eigenfunctions A s s we show are 
normalized to a common arbitrary value. In this Sub¬ 
section, we are not not considering the Coulomb contri¬ 
bution when we calculate the spectrum of Sp as a func¬ 
tion of temperature in Fig. 7. In the spin-degenerate case 
(J = OmHa) we see that the eigenvalues decrease in mag- 
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Figure 7. (color on¬ 
line) Spectrum of Sp(J = 
O.OmHa) as a function of 
T with only one negative 
eigenvalue that leads to a 
singular point. 



Figure 9. (color online) 
Eigenfunctions to a singu¬ 
lar eigenvalue at T « T cross . 
All eigenfunctions except 
one are of either type. 



t [ K ] 


Figure 8. (color online) 
Spectrum Sp(J = O.lmHa) 
as a function of T with 
many negative eigenvalues 
that cross zero and lead to 
singular points. 



Figure 10. (color on¬ 

line) A s s (t) at T c for 
J — O.OmHa, O.lmHa 
without the Coulomb 
interaction. 


nitude with temperature in a monotonous way. At low 
temperature all eigenvalues but one are positive valued; 
the negative eigenvalue crosses zero, at the temperature 
T C (J = OmHa) « 30K in the model present model, above 
which Sp becomes positive definite. 

As compared to the J = OmHa, the spectrum at finite 
splitting J = O.lmHa is fundamentally different. For 
small T we observe many negative eigenvalues and, most 
interestingly, several solutions det (Sp) = 0 at low tem¬ 
peratures. There is a temperature regime T cross « 10K 
(in this model) in which most negative eigenvalues cross 
zero and become positive. Beyond T cross only one negative 
eigenvalue remains, crossing later at T C (J = O.lmHa) ~ 
25K. Continuously reducing the splitting, this specific 
eigenvalue/eigenfunction pair can be traced to the spin- 
degenerate limit, although we do not show this here. Sim¬ 
ilarly, upon reducing the splitting, it is found that the 
temperature range where the crossings appear goes to 
7 C r OSS ~^ OK as J —> 0. 

We analyze the eigenfunctions in Fig. 9 and 10 corre¬ 
sponding to these multiple solutions and see that only the 
one at T « 25K has a continuous behavior. The other 
solutions are of two kinds and both show numerical dis¬ 
continuities (see Fig. 9). While one (green in Fig. 9) has a 
l/^-like pole, there is a second kind (red in Fig. 9) which 
has a delta peak like structure, i.e. the value at the pole 
of the first kind is large while the rest is extremely small. 


Increasing the number of sampling points increased the 
relative value at the discontinuity so this lead us to the 
conclusion that we are numerically sampling an unbound 
function. It has to be understood that an unbound func¬ 
tion cannot be the linearized solution to an originally 
non-linear fixed-point problem. This because at the pole 
the function is not small and a linearization cannot be 
justified. We expect that in the non-linear equation these 
type of solutions will be suppressed. We therefore ignore 
these other solutions in the following discussion and al¬ 
ways refer to the continuous, bound, high temperature 
solution. 

As a side remark we point out that comparing 
Tscdft^j _ 0 ) ~ 30K with the solution the Eliashberg 
equations, the latter predicts a much higher J , c Eliash (J = 
0) « 50K. For a detailed comparison, see Fig. 4 b) where 
we show the linearized i n the phase diagram of 

the Eliashberg equations. We also observe via the black 
curve of Fig. 4 b), that the Eliashberg solutions predict a 
SC phase that is less susceptible against a splitting. The 
reason for the lower T c prediction is that within the xc- 
potential construction G was replaced with G KS which 
violates Migdal’s theorem 11 . The solution has recently 
presented by Sanna et al. 25 using a corrected self-energy 
in the functional construction. We will come back to this 
point and elaborate on the distinction in the Appendix 
A. As a curious result, the linearized j , c SCDFT (j) curve 
bends upwards and starts an almost linear increase at 
the point where the transition is expected to become of 
discontinuous type. We investigate this issue in the next 
Subsection VA3. 


3. Analysis of the B 0 dependence of T c 

Using the condition of Sp to be positive definite we 
compute the T C (J) curve of the model (see Fig. 11). At a 
low field the T c ( J) curve behaves as expected; the critical 
temperature is slowly reducing with increasing J. Similar 
to the Eliashberg results in Sec. IV, the SpinSCDFT pair 
potential seems to be more resistant against a splitting 
than the BCS approach predicts. 

In the regime of a first order phase transition, where 
the conditions for a linearization are not met, Spin¬ 
SCDFT behaves differently as compared to the linear 
BCS solution of Fig. 3. While in neither case, BCS nor 
SpinSCDFT, a linearization can be expected to yield sen¬ 
sible results for a discontinuous first order transition, the 
behavior of the T C (J) curve from SpinSCDFT is certainly 
more unphysical. At high field, past J ~ 0.15mHa the 
T c ( J) curve bends outwards and starts an almost linearly 
increase with J. 

In figure lib), the eigenfunctions to the singular eigen¬ 
values of Sp for increasing J are plotted. We can clearly 
observe that the upturn the T C (J ) curve is accompanied 
by an increasing localization of Al at the Fermi level. The 
usual high energy tail gets more and more suppressed. 
For very large splittings, A s s becomes numerically noisy. 






















(a) T C (J) excluding the Coulomb 
coupling also removing triplet 
self-energy parts (singlet) or not (full) 
plus a BCS fit to the same T C (J = 0). 



(c) T C (J) including the Coulomb 
coupling together with a BCS fit to 
the same T C (J = 0). 


0 

^0 

<-0 

-0 



(e) T C (J) excluding the Coulomb 
coupling and Sp together with the 
linear BCS curve. 



(b) Normalized Z\g(e) a t E c full (J) °f 
(a). 


(d) Normalized Ag(e) a t E c full (J) °f 
(c). 


(f) The normalized A s s (z) from 
Sp - gf at T C (J) of (e). 


Figure 11. (color online) T C (J) from SpinSCDFT and the linear BCS curve with the same T C (J = 0) . We show the T C (J) in 
- or excluding the Coulomb coupling, in - or excluding contributions from triplet self-energy parts. Beyond the range of the 
continuous transition, the SpinSCDFT solutions start to increase at J ~ 0.15mHa. The effect of the Coulomb potential is to 
reduce the T C (J = 0) but the overall shape remains essentially unaltered. The normalized eigenfunctions are shown with a 
color code indicated the respective splitting. The solutions become numerically noisy at large splittings. 


Reintroducing the Coulomb coupling to Sp we observe 
a similar behavior. In this case A s s shows a character¬ 
istic negative tail induced by the Coulomb renormaliza¬ 
tion mechanism 26,27 as it occurs within SCDFT 11,12,28 . 
From the comparison between j , c SpinSCDFT (J) with the 
green dashed BCS curve in Fig. 11, we note that in the 
second order regime T BCS (J) scales down with J faster. 
In order to make the strong coupling SpinSCDFT the¬ 
ory more similar to the weak coupling BCS approach 
we disregard Sp in Fig. 11 e) and f). In this case we 
are only considering the effectively attractive coupling 
among electrons via phonons, similar to Frohlich 29 and 
BCS. The effective Frohlich interaction requires the cou¬ 
pling to be small, and moreover we neglect the phonon 
influence on the normal state (Nambu diagonal) part 
of the self-energy entirely. Thus this approximation is 
called the weak coupling limit. As expected, the re¬ 
sulting r SpinSCDFT ( J)/T SpinSCDFT (0) and T BCS (J)/T BCS (0) 
behave very similarly. Here the T C (J) curves shown in 
Fig. 11 e) also feature the linear increase for high split¬ 
ting. Moreover we observe a discontinuous jump of the 
critical temperature at a certain splitting J c which is ac¬ 
companied by the eigenfunction dramatically changing 


shape. After the jump, the solution does not have a 
common sign convention but shows positive and nega¬ 
tive parts. Also here we find numerically noisy solutions. 

The BCS T C (J) curve, fitted to the same T c (0mHa), 
matches the weak coupling SpinSCDFT T C (J ) curve 
Fig. 11 e), not the strong coupling curve of Fig. 11 a). 
This points out that the strong coupling Sp term does 
not simply scale T C (J) down equally on both, T and J 
axis. Instead, Sp leads to a larger T c ( J) reduction of the 
temperature axis. Thus we conclude that strong coupling 
systems are less effected by an exchange splitting relative 
to their T c (0mHa). 


B. Non-Linear Sham-Schliiter Equation 

The previous section has shown the importance to con¬ 
sider the fully nonlinear Sham-Schliiter equation Sp [Z\g] • 
zAg = 0 of Eq. (6) when working in the limit of strong 
external field/large exchange splitting J. We solve the 
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Figure 12. (color online) J — T diagram of solutions to the 
non-linear gap equation. We include the Tc Ull (J) curve 
(dashed blue) from the linearized functional of Fig. 11a). 


Figure 13. (color online) The SC gap in the SpinSCDFT 
Go Wo DOS. The dashed blue line is the linear Tc Ull (J) of 
Fig. 11 a). 



o.i 

£ [mHa] 




Figure 14. (color online) Z\g(e) f° r J — 
O.OmHa as a function of T. For low T, 
Z\g(e) goes to zero at e ~ 0. 


Figure 15. (color online) Z\g(e) f° r J — 
O.lmHa as a function of T. For low T, 
Z\s(e) remains above J at £ ~ 0. 


Figure 16. (color online) Z\g(e) f° r 
T — 10K as a function of J. At J ~ 
0.17mHa zlg(e) dramatically change 
shape. 


fully non-linear Sham-Schliiter equation 

Al = IC s [Al}-Al (9) 

JC s = S^-(S p + S), ( 10 ) 

with the splitting matrix S chosen to be S™(e, J = 
O.OmHa) (more details on this procedure can be found 
in Sec. I.C.2). In Fig. 14 we show results, neglect¬ 
ing the Coulomb coupling along the iso-splitting line 
J = O.OmHa as a function of temperature T. We ob¬ 
tain a Z\g(e) that goes to zero at the Fermi level for 
low temperatures (the purple to blue lines in Fig. 14). 
This means the SC KS system is not gapped (still main¬ 
taining x ^ 0) and we cannot directly interpret the SC 
KS excitations as quasi particles. In order to have the 
computationally convenient DFT scheme and a good ap¬ 
proximation to the quasi-particle structure at the same 
time we introduce the one-cycle Dyson equation itera¬ 
tion for SC in the Appendix A. This approach is similar 
to the common Go Wo approximation in band-structure 
theory 15 and leads to excellent results in SpinSCDFT. 

To complete the discussion of the J and T depen¬ 
dence of SpinSCDFT, we need a characteristic number 
of a given A a s (z) solution. As mentioned, Z\g(£ = 0? J) * s 
not a sensible choice, because it neither corresponds to 
an excitation gap nor is it a measure for the size of the 


potential Z\g(e). 

Instead, we chose f Z\g(e)d£ and the resulting Spin¬ 
SCDFT J — T diagram of Fig. 12 shows a transition at a 
point where, from the shape of the non-linear BCS and 
Eliashberg diagram the first order phase transition is to 
be expected. However, following this discontinuous tran¬ 
sition, the solutions Z\g(e) do not vanish but have a differ¬ 
ent shape. In Fig. 16, we show the Z\g(e) with increasing 
splitting on the equal-temperature line at T — 10K and 
the transition is clearly seen. In general, while before a 
critical splitting J C (T) the potential is little effected by 
the splitting, past J C (T) the solutions Z\g(e) localize at 
the Z\^(J) Fermi level and show positive as well as neg¬ 
ative regions. This behavior is similar to the shape of 
the potential from the linearized Sp as given in Fig. 11 
b). We show the T C (J) curve from the linear equation 
as a dashed blue line in Fig. 12 and see that it marks 
the border of the appearance of the curious solutions in 
the non-linear equation past the range in J of the second 
order phase transition. 

Due to the Coulomb renormalization, including the 
Coulomb repulsion, f Z\g(e)d£ is predominantly negative. 
Thus, as a physical property, we compute the number of 
condensed electrons N sc = f dr f dr'|%(r, r')\ 2 instead. 
We show the SpinSCDFT J — T diagram including the 
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J [mHa] 


Figure 17. (color online) N SC (T,J) and the linear T C (J) 
(dashed blue) including the Coulomb repulsion. The 
Z\®(e) after the transition lead to almost no condensed 
electrons. 



Figure 19. (color online) A s s (z) at T — 3K for several J 
along the dashed green line in Fig. 17. The Z\g(e) past 
the transition at J — 0.056mHa are similarly confined to 
the Fermi level region as without the Coulomb repulsion. 

Coulomb coupling in Fig. 17. The region past the tran¬ 
sition has essentially no condensed electrons, while, still, 
the A a s (z) i s n °t zero (compare Fig. 19). The SC Go Wo 
gap is shown in Fig. 18 and is very similar to the results 
without the Coulomb interaction of Fig. 13. Again, we 
find only small features past the transition (compare also 
Fig. 20). 


C. Numerical calculation of the DOS from the 

G 0 Wo GF 

We compute the GF according to the SC GoWo scheme 
derived in the Appendix A. In detail, we solve the 
Eq. (A38) using the Eqs. (A23) to (A26) together with 
the Eqs. (A2) to (A5) for the model and couplings de¬ 
scribed in Section II. We exclude the Coulomb potential 
at this point for a better comparison with Eliashberg the¬ 
ory although there is no difficulty to include it. 

In Fig. 13 we compute the GqWo corrected DOS at 



J [mHa] 


Figure 18. (color online) The SC gap in the SpinSCDFT 
Go Wo DOS. The dashed blue line is the linear Tc Ull (J) of 
Fig. 11 c). 



Figure 20. (color online) Up-spin channel of the GoWo 
DOS along the dashed green line in Fig. 18 at T = 3K. We 
see only small features from the unphysical, oscillatory 
solutions past the transition (light blue to red in Fig. 19.) 

every point in J and T and extract the SC excitation 
gap. We find that the curious solutions past the transi¬ 
tion J C (T) lead to almost no excitation gap. The reason 
is that for the self-energy in the calculation of the SC 
DOS in the Appendix A Z\g(e) * s integrated in e. If the 
high 6 region, away form the Fermi level are strongly sup¬ 
pressed, as in the KS potential past the J C (T), the effect 
on the excitation gap is negligible. 

Comparing with SpinSCDFT GoWo gap of Fig. 13 
with the BCS (Fig. 3) and the Eliashberg J — T dia¬ 
gram (4) we conclude that the point of the transition 
can be clearly identified. Moreover this one-cycle correc¬ 
tion sheds light onto the appearance of the Fermi-level 
localized solutions past the critical field J C (T). We have 
seen that for small T and J = 0 the non-linear zA®(e) 
go to zero at the Fermi level (compare Fig. 14) while 
the analogue of GoWo GF, the excitation gap of Fig. 13, 
takes its largest value at T = 0 and shows the expected 
monotonous decay with temperature to T c . 

This implies a significant difference in the quasi parti- 
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Figure 21. (color online) DOS from the Go Wo GF. We show the DOS in (a) corresponding to the SpinSCDFT results A®(e) 
shown in Fig. 14 for no splitting J — OmHa. In (b) (cr =^) and (c) (cr =t) we present the two different spin channels of the 
DOS with the SpinSCDFT results for the Ag(e) along the iso-temperature line T = 10K as shown in Fig. 16. 


cle states if a splitting occurs with such a 4\g(e). While 
the KS particle with the dispersion = sign(cr)J + 

ayje 2 + 4\|(e) 2 is strongly altered by the splitting be¬ 
cause the Bogoliubov branches change their order (com¬ 
pare the earlier discussion in the conclusion of Sec. V A 1 
and in I) this is not the case in the true quasi particle 
structure. In fact, from Fig. 15, we see that the SC so¬ 
lutions 4\®(e) if J > 0 do not go to zero and, instead, 
rise with J to prevent this situation. On the other hand, 
after the discontinuous transition we find 4\®(0, J ) < J. 

In the functional construction, the replacement G —>■ 
G KS is thus a strong suspect for the occurrence of this 
curious solutions past the SC transition. This is because 
G and G KS deviate in that the latter can be non-gapped 
while still corresponding to a SC solution. 


D. Triplet components 

The present implementation of SpinSCDFT assumes 
the spin decoupling approximation, i.e. assumes the pair¬ 
ing to be of spin singlet type (compare Sec. I.III.A.2.c). 
However, it was also shown in I that a magnetic splitting 
creates triplet components in the pairing potential, even 
for a purely singlet order parameter density. Triplet com¬ 
ponents appear as an intermediate step, in the self-energy 
that leads to the GO-functional in I since the Nambu off 
diagonal upspin and downspin components are in general 
not equal and of opposite sign. They can be intermediate 
since such triplet self-energy contributions lead to triplet 
as well as singlet order parameter contributions. The 
intermediate triplet self-energy that leads to singlet or¬ 
der parameter contributions can be included in the spin 
decoupling approximation functional without difficulties. 
From the theoretical side, this is an unpleasant signa¬ 
ture of formal inconsistency. We have, in fact, computed 
the critical temperature and KS gaps with and without 
these intermediate triplet self-energy terms. In Fig. V A 2 
and Fig. 11, we observe that their effect is negligibly 
small. The possibility of a triplet condensation, i.e. non¬ 


vanishing triplet order parameter contributions, in not 
investigated further in this work. 


E. Extension to real materials 

In this work, properties of the free electron gas with 
a phonon and Coulomb coupling subject to an homoge¬ 
neous exchange splitting have been calculated. To com¬ 
pute real materials without the use of adjustable param¬ 
eters, the electron-phonon coupling and the Coulomb po¬ 
tential has to be calculated from first principles. Then, 
according to the equations (1.122),(1.123) and (1.124) 
these couplings, as well as the computed single parti¬ 
cle states Ska may well have a distribution in J differ¬ 
ent from the homogeneous J 0 = —/i b B 0 that we are 
considering here. Also, sometimes, several regions in 
the Brillouin zone (or: in fc) have different couplings 
and a different SC pairing as in the well known case of 
MgB 2 30 . The isotropic formulation does not have to be 
given up, often it is enough to group this regions which 
we refer to as multi-band SC 30,31 . We extend notation 
e = (e:, J, b) f de = f ds f d J where b labels the 
groups of quantum numbers {k} sharing similar pairing. 


VI. SUMMARY AND CONCLUSION 

In this work, we have presented parameter free ab- 
initio calculations of a superconductor in presence of an 
homogeneous exchange splitting as for example the re¬ 
sult of an external magnetic field. We have used two 
approaches: A generalization of the Eliashberg approach 
and SpinSCDFT. SCDFT allows the direct inclusion of 
Coulomb interactions in a straightforward way, while its 
direct inclusion remains to be problematic within Eliash¬ 
berg where one has to rely on the /i* approach 32,33 . The 
Eliashberg equations, on the other hand, provide the ref¬ 
erence for the phononic self-energy, allowing to under¬ 
stand and develop functionals for SpinSCDFT. 
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We have implemented a code that solves the Spin- 
SCDFT equations with a linear and non-linear xc- 
potential and the non-linear Eliashberg equations derived 
in I. The xc functional is derived in I from the Sham- 
Schliiter equation based on the replacement of the in¬ 
teracting with the SC KS GF. We have investigated the 
behavior of the xc-potential on a model of a free electron 
gas with a tunable, homogeneous exchange splitting J, a 
phonon coupling that resembles to the one of MgB 2 and, 
optional, a static Coulomb interaction in the Thomas- 
Fermi approximation. We compute the SC properties of 
this system and find that in the regime of a second or¬ 
der phase transition in the T vs J diagram, SpinSCDFT 
results in a curve that compares similar in shape to the 
Eliashberg solutions. Removing the contributions in the 
functional that arise from the normal state (Nambu diag¬ 
onal) part of the self-energy we arrive at a shape that is 
very similar to the BCS behavior. Including the Coulomb 
interaction reduces the critical temperature but other¬ 
wise does not largely effect the shape of the J — T dia¬ 
gram. 

In agreement with BCS and Eliashberg, SpinSCDFT 
predicts a discontinuous transition in A s s (t) for large J 
except that the A a s (t) past the transition are not zero but 
have a curious shape that has positive and negative val¬ 
ues. Furthermore, the solutions Z\®(e) increasingly adopt 
non-vanishing values more or less only directly at the 
Fermi level £ ~ 0. In addition, we find that the non¬ 
linear SpinSCDFT solutions go to zero at 5 « 0 for T —>> 0 
and thus the SC KS GF is not gapped while the inter¬ 
acting and GoWo GF is. Since we have noted in I that 
the low center of energy range is where zA®(e) compared 
with J because the Bogoliubov eigenvalues at 5 = 0 read 
E%(0,J) = sign(cr)J + <a|Z\g(0, J)| we believe that this 
range is crucial. In contrast to the ones before the tran¬ 
sition, the curious A a s (z) P as t have 1^1 > |^g(0,«/)|. 

We perform a Go Wo like correction to the GF where 
we solve the Dyson equation with the same self-energy 
that we used originally for the xc-potential construction. 
The resulting excitation spectrum (here in the isotropic 
case the DOS) is gapped and behaves as one would ex¬ 
pect for a SC. From this result we conclude that a fitting 
technique of the self-energy similar to A. Sanna et a/. 25 
will allow us to reproduce the J—T diagram of Eliashberg 
while keeping the possibility to include the Coulomb po¬ 
tential in addition to a numerically simple form where the 
Matzubara summations can be computed analytically. 


Appendix A: Quasi-particle Excitations from the 
One-cycle Interaction Green’s Function 

The theoretical definition of SC is the existence of a 
non vanishing order parameter x (Eq. 1.4), while exper¬ 
imentally SC are usually characterized by the properties 
of their excitation spectrum, namely the single particle 
gap at the Fermi level 34 . This can be rather directly ex¬ 
tracted from the solution to the Eliashberg equations on 


the imaginary axis since A^ =0 (J) is closely related this 
excitation gap itself 24 and we use it in Fig. 5 for the J — T 
diagram. 

The SC KS system of SpinSCDFT is designed to re¬ 
produce the densities of the interacting system not the 
quasi particle spectrum. On the other hand, for a normal 
metal the KS particles are often in good agreement with 
experiment so that the resulting KS excitation spectrum 
is used as an approximation to the interacting quasi par¬ 
ticle spectrum. With the potential A a s (z) °f Eig. 14 it 
turns out in SCDFT, also for the zero field case 25 , this 
is not always the case, since e.g. for T —0 the SC KS 
system is not gapped. 

To predict a proper excitation spectrum without hav¬ 
ing to solve the Many-Body problem self-consistently we 
introduce the GoWo approximation in the context of SC. 
This means to solve the Dyson equation once while re¬ 
placing the interacting GF with the SC KS GF in the 
self-energy. Here we use the same approximations for the 
self energy made to arrive at the functional in I which 
means we use XJ KS = XJ[G KS ] instead of the true self- 
energy S[G]. 

In this Section we work in the isotropic formulation 
but note that the approach is easily generalized to the 
anisotropic case. We use the notation e = 6, J, b and the 
isotropic Dyson equation 

G„(e) = ((G-(e))- 1 + r-(e)) _1 1 (Al) 


that follows form the assumption that the couplings de¬ 
pend on k via the center of energy £fc t~£-fc,r £ anc [ ^h<e 

splitting £fc t~g-fc,r j anc [ e isotropic bands b (that 
is a set of quantum numbers {&}). We introduce the 
notation (e) = h*(t)G%£± k ± a (un)- The averag¬ 
ing procedure //^(e) on equal splitting and equal center 
of energy surfaces is defined in Eq. 1.120 . We refer to 
the non-vanishing matrix elements with a spin label that 
refers to the first index of G k ’J ±k±cr (uj n ) and similar for 
the self-energy. 


1. Imaginary Axis Formulation 

The inversion of the Dyson Eq. (Al) to compute 
the GF explicitly is very analogous to the derivation 
of the Eliashberg equations in I, Sec. I.IV.A.l . We 
compute G n (t) via Eq. (Al) and the non-vanishing 
components are found to be (suppressing the ar- 
guments of ff n<7 (e), i£(e), itf(e), Ag («), K* M, 
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and r;±(e)) 

q 1<1= 1 v- ff nCT + a(e+£^+sign((T)A)f) 

2 ?n<j “ iw n ——sign(u)(j + ^n) —tt$ncr 

(A2) 

g_!,_ 1= 1 ^n-a+aje+^r- sign(cr)A^) 

23n-<r ^ iw n -Z , “+sign(cr)(j+i;^)+aS r n _ <T 

(A3) 

- x ,^signOr) v “(^ A +i^ A +sign(a) (Z* n ~ + i£+)) 
Gn<T 2^ E iw n -r“-sign((7)(J+r^)-gJ raCT 

(A4) 

g _ ia _ sign(q) « ~ i^g A +sign(a) (Xfc- - i% + )) 

25?n-<r “ iw n -Z , “+sign(cr)(j+A , ^)+gS r „_ (T 

(A5) 


with 


to (1.69) are given by 


^ = fa fac?F^(t,t',n) x 

x E^l^W^£“>n) ( A15 ) 

a 

KCn' 1 = fa J d t'c?F°(ty,!2) X 

( A16 ) 

a "" 

= —sign(cr) fa fa c?F(t, t', 12) x 
*Y,°§r M ^ E °^) (Ai7) 

a 

KCn 1 = - si S n (^) fa fac?F(e,e',f2) X 

rv/\ s/ * 

x E^r M ^(^- jB “>-) ( A18 ) 


S'na(e) = ((e + + sign(cr)A^) 2 + 

+ (^ A + i^ A + sign(cr)(i;; + + r;-)) x 

x (rr - i^r + sign(<7)(i;;+ - r;-))) 5 

(A6) 

where the self-energy parts are constructed similar to the 
Eliashberg theory with the result 

^^E^+^r 1 ’ -1 ) (A 7 ) 

a 

AT «1 E si S n ( <T ) (^n 1 ’ 1 + 1 ’" 1 ) (A8) 

<J 

Z £ n = fa{Kn 1 ’ 1 -Zan- 1 ’- 1 ) ( A 9) 

a 

= \ E si g n ( CT ) (VSl 1 ’ 1 - ^n -1 ’” 1 ) ( A1 °) 

a 

k + = \ E (^ 1,_1 + ^r 1,1 ) ( An ) 

<j 

K~ = \ E(^™ 1 ’“ 1 - A - _1,1 ) ( A12 ) 

a 

+^r 1,1 ) ( Al3 ) 

a 

= \ - ^r 1 ’ 1 ) • (A 14 ) 

<7 

Note, however, that E^ 1,_1 contains a triplet contribu¬ 
tion that is generated by the coupling imbalance of the 
spin channels. The isotropic variants of the Eqs. (1.66) 


with A S J short hand for Z\g(e'), the averaged A s sk and 

F' = yje' 2 + A s s ' 2 . Furthermore Efa = sign(cr)J + 
sign(<a)F' and similarly the Eqs. (1.75) and (1.76) be¬ 
come 

i^^-sign^E fde'^C“‘“(e, z')fp{ET) 

(A19) 

£cTn=sign(a)£ /de'^-C““‘(e, t')fp{-E?). 

a ^ 

(A20) 

With these equations we can compute the Go Wo GF from 
the results of a converged SpinSCDFT calculation. 


2. Real Axis Formulation 

To obtain the (L)DOS from the temperature GF we 
substitute 

icj n oj + irj (A21) 

where r] is a real positive infinitesimal 35 . The expres¬ 
sion Eqs. (A2) to (A5) remain essentially unchanged on 
the real axis, except that we have to insert the SE parts 
Eq. (A7) to (A14) on the real axis and write i 77 +cj instead 
of the Matsubara frequency. Here we have two options, 
first we may compute the SE parts on the imaginary axis 
and use a numerical analytic continuation to the real axis, 
or we can compute analytic formulas for the real axis and 
use them. We choose the latter because this avoids the 
sometimes unstable analytical continuation. 

We will see that the SE parts, e.g. E^ A (e), on the real 
axis have to be computed via independent calculations 

















14 


of imaginary and real part. The dependence on the Mat- tion in the SE. Thus on the real axis 

subara index of the SE is only via the function M ph of (O) f (F) (n\ + f (—j?\ 

Eq. (1.73), i.e. the results of the first Matsubara summa- M ph (J2, E,uS) = P ^ —- — P ^ — — — -- 

-m((np(ft) + fp(E))5(ft -E + ui) 

+ (np(ft)+fp{-E))b(ft + E-uo)). (A22) 


Here P is the principle value operator. Because of the 
very different nature of the imaginary and real part of 
the SE we compute both parts independently. Then we 
obtain 


3£“(e, w) = -7T Jdt 1 J2 ^W~((M E ? - w) + ME?)) (a?J£(e, t',E? - w) - e> - E?)) + 


(ME? + oo) + fp(Ep) (a 2 F M D (e, e', E? + w) - c?F°(c, e', -E? - oo 


^{t,oo) = Jdft 


r MM) + fp{E?) 
Q — E^f oo 


fia A 6 

■M n ) + MO fiM® + M~E?) , ~np(ft) + fp{-E?) 


ft - E? - oo 


ft + E? - oo 


+ P- 


ft + E? + oo 


) 


3Z" 


\c,co) = -n jdt' E ((ME? - w) + ME?)) (a^e, e', E? - oo) - a 2 F^(t, e> - E?)) 

J fia 

-(ME? +co) + ME?)) (cPF^e, e', E? + oo)~ c?F°(c, e', -E? - oo 


^(z,oo) = Jdft fdc'J2^jyEa 2 F°(z,z',n)(] 


t MM + ME?) 
ft - E? + 00 


+P 


Ha A 6 

np(ft) + fp(E?) ~np(F!) + M-E?) ~np{ft) + fp{-E?Y 


ft - E? - oo 

h 1 


-P- 


f2 T E?f — oo 


-P- 


Q + E^ f + oo 


(A23) 


(A24) 


(A25) 


(A26) 


and very similar for A OJZ (too ) that only differs from ZJ U 
by putting a sign(/i) into the spin sums. We also obtain 
E J (too) from the relation for E s (t,oo) in the same way, 
i.e. we put a sign (fi) into the spin sum. The above equa¬ 
tion again points out the problem in the s' integral if 
the energy dependence of d 2 F / ^(^e / , i?) is neglected. Here 
E™' —» a\s'\ for large \s'\ so there are parts in the inte¬ 
gral that behave as ^ leading to logarithmic divergence. 
Thus we see explicitly that we cannot compute the energy 
renormalization without considering the influence of the 
interaction on the full energy spectrum and quasi-particle 
occupations as was already discussed in I and Ref. 11. 


We define the integrand 


9S±(e, t', oo) = 

TT^sign ^ s MyAl( M E?-oo) + ME?)) X 

fia 

x (a 2 F M D (e, e', E? - oo)-c?F^(z, e> - E?)) (A27) 

= - E ( kfta 2 F-(z, e', ft) x 

f np(ft) + ME?) ^np(ft) + fp(—E?)\ 

X \ f2-E?+uo ft + E?-oo ) + 

+ME?)C^(t,t’)) 


(A28) 
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and further introducing 

B‘ n (t) = 2™(z)+ irr(e) (A29) 

K(e) = (e) - i^ A (e) (A30) 

B*(e) = I£ + (e) + I£-(e) (A31) 

B‘*(e)^r; + (e)-r;-(e) (A32) 

we obtain the following equations on the real axis 

B s (e, u)=J d t'A'jB- (e, e', w) (A33) 

B s *(e,w) =Jdt'AfB-(e,e',u]) (A34) 

B\z,uj)Mdt'A a s 'B + (t,t',uj) (A35) 

B t *(e, w) =Jde'AfB+(e,e',u) (A36) 

and thus, Eq. (A6) becomes on the real axis (omitting 


the arguments z,u) 

da = ((e + £ e + sign (a)A uz ) 2 

1 

+ (B S + sign((j)i? t ) ( B s * + sign(<r)B**)j 2 .(A37) 

Now we can finally obtain the retarded GF with the equa¬ 
tions from Eqs. (A2) to (A5) together with Eq. (A6) for 
fi neT (e) in terms of B and the corresponding SE parts con¬ 
structed from real and imaginary part close to the real 
axis. Then we can evaluate the DOS according to 

Paa{u) = -2 /de3(lim lim G“^(c))e(c) (A38) 

We obtain the local DOS p aa (r,oo) simply by re¬ 
placing g(e) with the local double DOS p a (e,r) = 
ika(e)\ifi k (ra)\ 2 . 
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